The routine solves a general Gauss-Markov linear model (GLM) problem:

where A is an N-by-M matrix, B is an N-by-P matrix, and d is a given N-vector. It is assumed that M <= N <= M+P, and 

Under these assumptions, the constrained equation is always consistent, and there is a unique solution X and a minimal 2-norm solution y, which is obtained using a generalized QR factorization of the matrices (A, B) given by 

In particular, if matrix B is square nonsingular, then the problem GLM is equivalent to the following weighted linear least squares problem

where inv(B) denotes the inverse of B. The sign _2, denotes Norm L2. 

References: 1.) Lapack v3.4 source code

Matrix Methods